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Abstract 

We study the first passage time (FPT) problem in Levy type of anomalous 
diffusion. Using the recently formulated fractional Fokker-Planck equation, 
we obtain an analytic expression for the FPT distribution which, in the large 
passage time limit, is characterized by a universal power law. Contrasting 
this power law with the asymptotic FPT distribution from another type of 
anomalous diffusion exemplified by the fractional Brownian motion, we show 
that the two types of anomalous diffusions give rise to two distinct scaling 
behavior. 
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I. INTRODUCTION 



For a stochastic process, the first passage time (FPT) is defined as the time T when 
the process, starting from a given point, reaches a predetermined level for the first time, 
and is a random variable [^. Escape times from a random potential, intervals between 
neural spikes, and fatigue failure times of engineering structures are all examples of FPTs, 
arising in physics 0, biology [Q, and engineering respectively. Thus, knowledge of the 
FPT distribution, f{t), is essential for the effective application of probabilistic analysis. (As 
a convention we use capital letters to denote random variables and lower case letters to 
denote their values.) Unfortunately, only in very few cases does one have explicit analytical 
expressions for f{t). One such case is the ordinary Brownian motion, an example of ordinary 
diffusion, in which the FPT is described by the famous inverse Gaussian law The main 
contribution of this work is the derivation of the exact solution of f{t) for a much broader 
class of stochastic processes, namely, the Levy type of anomalous diffusion in which the 
mean square displacement of the diffusive variable X{t) scales with time as < X'^{t) >r^ V 
with < 7 < 2. Specifically, using a recently formulated framework of fractional Fokker- 
Planck equation (FFPE) we express f{t) in terms of Fox or H- functions [§,0, which 
is shown to contain the inverse Gaussian distribution as a special case. Furthermore, we 
show that in the asymptotic limit of large t, /(t) scales with t as f{t) ~ t^^^'^/^. Our next 
result concerns the comparison with a different type of anomalous diffusion, represented by 
fractional Brownian motion (fBm) [^,|T^, where again < X^(t) >~ V with < 7 < 2. For 
this type, we argue that f{t) ~ t^l'^^'^ for large t, a result that has been conjectured earlier 
|r3| . Finally, we present numerical simulations which verify the analytical results. 



II. FFPE AND DERIVATION OF FPT DISTRIBUTION 

The Levy type of anomalous diffusion considered in this work is a class of non Gaussian 
and non Markovian processes founded on the continuous time random walk (CTRW) where 
the waiting time obeys certain power law distribution 0. Let (j){y,u) denote the joint 
probability density between the waiting time U and the jump size Y. It can be shown 
that, depending on the specific form of (j){y,u), the CTRW can produce both subdiffusive 
(0 < 7 < 1) and superdiffusive processes (1 < 7 < 2) as well as ordinary diffusion (7 = 1) 
PJT^. For example, consider 

(j){y, u) = -jL= exp[-yV2a2li^^^4^, (1) 

where Y and U are decoupled with Y being a Gaussian variable. (We note that, strictly 
speaking, the distribution of U is not a Levy stable distribution, but belongs to the domain 



of attraction of a one-sided stable Levy law. We call U a "Levy type of variable" for 
want of a better name.) For 1 < a < 2, the corresponding CTRW is characterized by a 
subdiffusive process with 7 = a — 1, and for a > 2, one gets ordinary diffusion with 7 = 1. 
If, on the other hand, Y and U are coupled through 

<p{y, u) = UiulT - |y|/a)-^^-=i^, (2) 
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where 2 < P < 3 and S{-) is the Dirac delta function, the CTRW describes a superdiffusive 
process with 7 = 4 — /?. 

Let W{x,t) be the probabihty density function for a CTRW X{t) with X(0) = 0. Con- 
sider the generahzed diffusion hmit where a and r are scahng parameters for the space and 
time variables. For the subdiffusive case, this means taking the limit and r ^ 

with K = o"^/2r(l — 7)r'^ kept a constant, and for the superdiffusive process, this means 
taking the same limit with K = (3 — 7) (2 — 7)r(7 — 1)(T^/2(5 — 2'j)t"' kept a constant. In 
this limit it can be shown that the evolution of W{x,t) is determined by the following 
FFPE §: 

Wix,t)-W{x,0)= ,DT^K—W{x,t), 0<7<2, (3) 
where the Riemann-Liouville fractional integral operator qD^'^ is defined as [|TR|TB 

oD;^W{x,t) = ^J' dt' {t - tr-'W{x,t'), 7 > 0, (4) 

r(7) Jo 



with T(z) being the gamma function |Tj]. The constant K is the generalized diffusion 
constant defined in the above generalized diffusion limit. From Eq. (|^) it is easily shown 
that < X^{t) >= 2Kt^/T{l + 7). 

In the framework of the Fokker-Planck equation, the first passage time problem is defined 
in terms of having absorbing boundaries at a; = —00 and x = a, where a is the predetermined 
level of crossing, with the initial condition W{x,0) = 6{x) An equivalent formulation, 
due to symmetry, is to solve Eq. with the following boundary and initial conditions: 
W{0,t) = 0, W{oo,t) = 0, W{x,0) = 6{x — a), where a; = a is the new starting point of 
the process, containing the initial concentration of the distribution. (This latter formulation 
makes the subsequent derivation less cumbersome.) Once we solve for W{x,t), the first 
passage time distribution f(T) is given by |1| 

f(i) = -Jtl dxWix,t). (5) 

Taking into account of the boundary and initial conditions we are led to the following 
expansion for W{x,t) |ir 
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W{x,t) = — dksmkx sinka A{k,t), (6) 

71 Jo 

with A{k,0) = 1. To determine the unknown function A{k,t), we substitute the above 
expansion for W{x,t) in Eq. and, after straightforward algebra, obtain A{k,t) — 1 = 
—Kk'^ QD^^A{k,t). Taking the Laplace transform with respect to t, we have 



where A{k^p) is the Laplace transform of A{k,t). Here we have applied the result |]TB[ that 
the Laplace transform of qD^^ A{k^t) is A{k,p)/p'^. Inverse Laplace transform of Eq. 
yields PI 
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A{k,t) = E^{-k'^Kt'^), 



(8) 



where Er^^z) is the Mittag-Leffler function [T^. Substituting Eq. (||) into Eq. (||) we get 

2 /""^ 

W{x,t) = — / dksinkx sinka E^^—k'^Kf). 
n Jo 



(9) 



To proceed further, we introduce the Fox or H-function PJIOj] which lias the following 
alternating power series expansion: 



TTm,n 



i=l A:=0 



A;!R 



X 



UT=lJ^l ^jb, - B.Sik) Ur=l -ar + ArSi,) 

n«=m+l — + BuSlk) ]Xv=n+l r(ai. — 



(10) 



where s^fc = (&/ + k)/Bi and an empty product is interpreted as unity. Further, m, n, p, q are 
nonnegative integers such that < n < p, 1 < m < q; Aj, Bj are positive numbers; aj, bj 



can be complex numbers. The H-function has several remarkable properties |]T0| which are 
listed in the Appendix. 

By comparing the series expansion |jl9[ of the Mittag-Leffler function E^{z) with that of 
the H-function [cf. Eq. ([To|) ], Eq. (||) can be rewritten as 



W{x,t) 



2 1 1 f 

J dksinkx sin/ca Hi 2 ik Kf 



(0,1) 

(0,1), (0,7) 



(11) 



Letting k' = k{Kt^f/'^ and using Property 5 [Eq. (|^)] of H-functions, the above equation 
becomes 



W{x,t) 



27T{Kr') 



— —777 / dk'\cosk'(x — a) — cosk'(x + a)]H}'}, ik' 

:pyi^ Jo L V ^ \ ji 1,2 y 



(0,1/2) 

(0,1/2), (0,7/2) 



(12) 



The Fourier cosine transforms can be solved by successive applications of a Laplace and an 
inverse Laplace transform (a technique pioneered by Fox for solving a wide variety of 
integral transforms) to give pT| 



W{x,t) 



» I. 

Q 



\x — a\ 



1 rr2,l ( 3: + a 

2{x + a) ^^^\{Kt^Yl^ 



1,1/2), (1,7/2), (l,l/2)\ 
(1,1), (1,1/2), (1,1/2); 

:i,l/2), (1,7/2), (1,1/2) 
(1,1), (1,1/2), (1,1/2) 
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Now, applying Properties 2, 1, 3 and 6 of the H-functions (listed in the Appendix) in the 
given order, we finally get 



W{x,t) = 
1 



rrl.O 

-"1,1 



\x — a\ 



(1-7/2,7/2) 

(0,1) 



1,0 



X + a 

JmnyT^ 



(1-7/2,7/2) 

(0,1) 
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Substituting Eq. (0) into Eq. 
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fit) 



d 
dt 
d_ 



we have 



dx H\'i 



2(irn)V2 

1 I 



2{Kt 



\x — a\ 
X + a 



[1 -7/2, 7/2) \ 
(0,1) ) 

:i-7/2,7/2)\ 
(0,1) I 



(15) 



Removing the exphcit time dependence in the integrands by rewriting the integrals in terms 



oi z = (x - a)/{KT^) 
evaluated to give: 

fit) 



'7U/2 



(X 



+ a) / {KT'^y^'^ , the above integrals can be explicitly 
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7t1,0 



(1-7/2,7/2) 

(0,1) 



(16) 



which is the main result of this paper. It should be noted that H-functions were first used 
in the context of probability distributions by Schneider ||2^. They have also been used to 
express solutions of fractional diffusion equations [^. In addition, the FPT problem in the 



(17) 



context of Levy processes has been considered in P4| . 

The series expansion of the H-f unction in Eq. ( |T6|) [cf. Eq. (jl^)] is 



fit) 



07 



E 



2irV2t(2+7)/2 ^ k\r{l - 7/2 - ^7/2)' 



This turns out to be also the series expansion of the Maitland's generalized hypergeometric 
function or the Wright function Qipi [|l^. Thus, an alternative expression for f{t) is 



fit) 



07 



2i^l/2^(2+7)/2 



oi'i 



(1-7/2,-7/2)' (irT7)i/2 



(18) 



For 7 = 1 (ordinary Brownian motion), the Wright function reduces to the following simple 
formula 



fit) 



{AnKt^y/^ 



-a?/AKt 



(19) 



This is the expected inverse Gaussian distribution for the FPT distribution of the ordinary 
Brownian motion 0. 

Next, we consider the asymptotic behavior of the FPT distribution for large values of 



t. Refer to Eq. (|T6D. Let z = a/{KVy/'^. It is known that [|TO|,0, for small z, hI'^' 
|^|6i/-Bi _ since 61 = and Bi = 1. Therefore, the FPT distribution /(t), for large t, is 
characterized by the power law relation 

f{t) ~ t-'-^l\ (20) 



which becomes the well known —3/2 scaling law for the ordinary Brownian motion. This 
power law behavior has been observed earlier by Balakrishnan for subdiffusive processes 
(0 < 7 < 1) using a different method. Using our method the same scaling law is shown to be 
applicable also to superdiffusive processes. After some manipulation, we can also determine 
the location trnax of the maximum of the FPT distribution: 
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where d = {1 - 7/2)(7/2)^/(2-7)(^q//^)2/{2-7)^ Yiom Eq. we see that the mean first 
passage time and all higher moments of the FPT distribution are undefined for < 7 < 2. 

The theoretical prediction for the full FPT distribution given in Eq. ([T6|) is verified 
by numerically simulating the underlying CTRW process characterized by the probability 
density function (f){y,u) [cf. Eq. (Q)]. For the sake of numerical efficiency, we replace the 
waiting time distribution in (/>(?/, u) by the Pareto distribution which is well justified 
for small values of r. Ten million realizations of the CTRW process are used to generate 
the numerical FPT distribution. The results are shown in Figure 1 for 7 = 0.5, a = 1.0, 
r = 10~^ and K = 0.1. We note that the numerical simulation is in excellent agreement 
with the theoretical prediction. The agreement would get even better as the generalized 
diffusion limit is approached (that is, as r — > and cx^ — >■ with K held a constant). 



III. ANOMALOUS DIFFUSION OF THE FBM TYPE 



Fractional Brownian motion X{t) [|n|,|T^ is a Gaussian process with X{0) = 0, < 
X{t) >= and < [X{t) - X{s)f >= \t - (0 < 7 < 2). By definition it provides 
us with another type of anomalous diffusion. The exact FPT distribution of this process 
is not known. It was conjectured [0, based on scaling argument and numerical evidence. 



that for large t, f(t) scales with t as 



f{t) ~ t^/2-2_ ^22) 



Notice that this power law behavior is different from that in Eq. ( pOD even though the mean 
square displacement < X(ty > has the same power law behavior (< X(t)'^ >~ t^) for both 
types of anomalous diffusion. Below we give a heuristic argument for this power law using a 
recent result ||2^ concerning the distribution of the maximum of a fBm over a given interval. 



Without loss of generality we set the threshold at a = 1. Let the probability that the 
maximum Mt of the fBm X{t) with X{0) = is less than 1 in the time interval [0,t\ be 
denoted by Pit): 

P{t) = Prob(Mt < 1) (23) 

Clearly, P{t) is also the probability that the first passage time T of the fBm is greater than 
t: 

Prob(r >t) = Pit). (24) 
This implies that the first passage time distribution fit) is given by 

fit) = |Prob(T <t) = -|Prob(T > t) = -^P(t). (25) 



Recent work by Molchan |2^ shows that, in the large t limit, P(t) scales with t as 

Pit)r^t^/^-\ t^oo. (26) 
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Substituting this in Eq. (^) we obtain Eq. p^ ) 
We present numerical results to verify Eq. 



29|| for the power 



j). Sinai's formula 

spectrum of the fractional Gaussian noise (fGn) is used to generate the fBm. The log-log 
plot of the FPT distribution is shown in Figure 2 (7 = 1.5). It is clear that the predicted 
slope of 7/2 — 2 = —1.25 is in excellent agreement with the numerical simulation. 



IV. CONCLUSIONS 

The two types of anomalous diffusions considered in this work lead to two distinct scaling 



behavior, Eq. (pOD and Eq. (p2D, for the respective FPT distributions in the asymptotic 
limit, despite the fact that they are both described by the same mean square displacement. 
Equation (^) is expected to be applicable to all CTRW types of processes, regardless of the 
specific forms of (f){y,u), for which the generalized diffusion limit leads to Eq. @. On the 
other hand, we expect Eq. ( [2^ ) to hold for Gaussian processes where < [X(t) — X{s)]'^ >~ 
\t — s\'^ for large \t — s\. 

In this work we considered only processes with < X(t) >= where the asymptotic 
limits of the FPT distributions are described by power laws. For a fBm, little is known 
about its FPT distribution when < X(t) >^ 0. For a Levy type of diffusion process, 
some exact results can be derived for the Laplace transform of the FPT distribution when 
< Xif) >7^ 0. We will present these results in other publications. 
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APPENDIX: PROPERTIES OF H-FUNCTIONS 



The H-function has the following properties [JIUI which have used in the main text. 
Property 1. The H-function is symmetric in the pairs (oi, Ai), . . . , (a„, A„), likewise 



, yap, Ap)] in {pi, Bi), . . . , (6^, B^^) and in {bm+i, ^m+i), • • • , Bq) 



Property 2. Provided n > 1 and q > m 



(ai,Ai), (02,^2), 
(&i,5i), 

(02,^2), ■ ■ - , (ap,v4p) \ 
(61, (6,-1, )• 



Property 3. Provided m > 2 and p > n 



Trm-l,n 
-"p-l,g-l I ^ 



(ai,Ai), (ap_i,y4p_i), 
(61, Si), (62,52), (6„5,) 

(ai,Ai), ■ ■ ■ , (ap_i, Ap_i) \ 
(62,52), {b,,Bg) ■ 



Property 4- 



HI 



{bj,Bj)j=i^ 



TTn,m 



Property 5. For A; > 0, 



(flj, Aj)j=l^,,,^p 



Property 6. 



(Oj, Aj)j=l,...,p 



(&„5, 



J, ^j)j = l,...,q 



p,g 



(1 - 6j,5j)j = l„„,g 

(1 — a^, v4j)j=i^ p 



(oj, A;Aj)j=i^...^p 



(Oj + pAj, Ajr)j = l,...,p 
(6j + pBj,Bj)j = i^,„^q 



(Al) 



(A2) 



(A3) 



(A4) 



(A5) 
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APPENDIX: FIGURE LEGENDS 



Figure 1: Comparison of tlie tlieoretical FPT distribution (solid line) with the distribution 
(dashed line) obtained by numerically simulating the underlying CTRW process for a 
Levy type anomalous diffusion with 7 = 0.5. 

Figure 2: Comparison of the theoretically predicted power law behavior of the FPT 
distribution for a fBm with 7 = 1.5 (dashed line) with the numerical simulation (solid 
hne) . 
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